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ABSTRACT 

Line-driven winds from stars and accretion disks are accelerated by scattering in numerous line transitions. 
The wind is believed to adopt a unique critical solution, out of the infinite variety of shallow and steep solutions. 
We study the inherent dynamics of the transition towards the critical wind. A new runaway wind mechanism is 
analyzed in terms of radiative-acoustic (Abbott) waves which are responsible for shaping the wind velocity law and 
fixing the mass loss. Three different flow types result, depending on the location of perturbations. First, if the 
shallow solution is perturbed sufficiently far downstream, a single critical point forms in the fiow, which is a barrier 
for Abbott waves, and the solution tends to the critical one. Second, if the shallow solution is perturbed upstream 
from this critical point, mass overloading results, and the critical point is shifted inwards. This wind exhibits a 
broad, stationary region of decelerating flow and its velocity law has kinks. Third, for perturbations even further 
upstream, the overloaded wind becomes time-dependent, and develops shocks and dense shells. 
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1. Introduction 

Radiation-driven winds which are accelerated by ab- 
sorption and rc-cmission of continmini photons in spec- 
tral lines form an interesting class of hydro dynamic 
flows, termed line-driven winds (LDWs). They occur in 
OB and W-R stars, cataclysmic variables, and proba- 
bly in active galactic nuclei and luminous young stellar 
objects. These winds are characterized by a unique de- 
pendence of line force on the velocity gradient in the 
flow. This causes a new, radiative wave type. 

The nature of these waves was first discussed by Ab- 
bott (1980, hereafter 'Abbott waves'), who found two 
modes, a slow, acoustic one propagating downstream, 
and a fast, radiative one propagating upstream. The 
critical point found by Castor, Abbott, & Klein (1975, 
hereafter CAK) from analysis of the stationary Euler 
equation is a barrier for these waves, the same way 
as the sonic point is to sound waves. The intriguing 
property of Abbott waves is that they propagate down- 
stream slower than the sound speed, while they prop- 
agate upstream at very large speeds, highly superson- 
ically. The latter fact reflects essentially the radiative 
nature of Abbott waves. 

No such upstream propagating, radiative mode was 
found by Owocki & Rybicki (1986), who calculated the 
Green's function for winds driven by pure line absorp- 
tion. The explanation is that the pure absorption case 
suppresses the radiative upstream mode, as photons 
can propagate only downstream. The radiative up- 
stream mode is back when line scattering is included 
(Owocki & Puis 1999). 

The question arises for the physical interpretation 
of Abbott waves. Do they represent a physical entity 
which is responsible for shaping the flow by communi- 
cating essential flow properties between different points 
in the wind? In particular, could Abbott waves be the 
prime cause for evolution of LDWs towards a CAK- 
type, steady-state solution? 

The analysis by CAK of the steady-state Euler equa- 
tion for LDWs has revealed an infinite family of math- 
ematical solutions, but only one, hereafter 'critical so- 
lution,' which extends from the photosphere to arbi- 
trary large radii. Other solutions either do not reach 
infinity or the photosphere. The former solutions are 
called shallow and the latter ones — steep. The unique, 
critical wind starts as the fastest shallow solution and 
switches smoothly to the slowest steep solution at the 
critical point. 

The shallow wind solutions found by CAK arc the 
analog to solar wind breezes, in that they are sub- 
abbottic everywhere. They were abandoned by CAK 
because they cannot provide the required spherical ex- 
pansion work at large radii. This exclusion of shal- 



low solutions can be criticized in different respects: (i) 
the breakdown happens only around 300 stellar radii, 
where basic assumptions of the model (fluid descrip- 
tion; spherical symmetry; isothermality) may become 
invalid, (ii) Shallow solutions could be extended to 
inflnity by jumping at some large radius to a deceler- 
ating wind branch. The latter was excluded a priori by 
CAK. However, this jump can occur beyond a few stel- 
lar radii, where the wind has already reached its local 
escape speed, (iii) For models of disk LDW, even the 
critical solution itself does not extend to infinity, and 
becomes imaginary beyond a certain radius (Feldmeier 
& Shlosman 1999). Jumps to the decelerating branch 
are unavoidable then. 

With shallow solutions being valid stationary solu- 
tions, the question arises, what forces the wind to adopt 
the critical CAK solution? Numerical aspects of this 
question have been discussed by Feldmeier, Shlosman 
& Hamann (2001), who noted that outer boundary con- 
ditions and a Courant time step which do not account 
for Abbott waves can set off numerical runaway, often 
towards the critical solution. 

In the present paper, we focus on a physical in- 
terpretation of Abbott waves, extending our previous 
work on this subject (Feldmeier & Shlosman 2000). We 
show that these waves are the prime driver of evolu- 
tion of LDWs towards a unique, steady-state solution 
characterized by a specific velocity law and mass loss 
rate. In particular, we find that as is the case for so- 
lar wind breezes, shallow solutions can evolve due to 
waves which propagate upstream to the wind base (the 
photosphere). As a new effect in LDWs, Abbott waves 
'drag' the solution in one preferred direction, towards 
larger velocities. The wind becomes stable when a criti- 
cal point forms, through which outer perturbations can 
no longer penetrate inwards. 

2. Abbott waves 
2.1. Wind model 

Only wind acceleration due to a line force in Sobolev 
approximation for radiative transfer is considered in 
this paper. The large number of lines driving the wind 
is dealt with using a CAK line distribution function. 
The latter is characterized by a power law index a, 
which lies between and 1. The Sobolev force is pro- 
portional to 

51 ~ y" dwn/„r^", (1) 

where n is the unit vector pointing in the direction of 
the surface angle element div, I is the frequency inte- 
grated specific intensity, and r is the Sobolev optical 
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depth, 



n.(n.(Vv))' 



(2) 



with rate-of-strain tensor Vv, mass absorption coeffi- 
cient Av, and ion thermal speed ?;th- The force propor- 
tionahty constant is fixed later in terms of the critical 
solution. Through Vv, the line force depends on the 
streamline geometry. To simplify, we take t in (|l|) out 
of the integral and replace it by an average or equivalent 
optical depth in one direction. Specifically, the latter is 
assumed to be the flow direction. This approach corre- 
sponds to the CAK 'radial streaming approximation.' 
For a planar flow with height coordinate z, the line 
force becomes, assuming constant k and wth, 



51 ^ F{v'/pr, 



(3) 



where v' ~ dv/dz, and the radiative flux F is a function 
of z. Even this highly idealized line force depends in a 
non-linear way on the hydrodynamic variables v' and 
P- 

So far, the following assumptions were introduced: 
(i) the line force is calculated in Sobolev approximation, 
using (ii) a CAK line distribution function (without 
line overlap), and (iii) applying the radial streaming 
approximation for (iv) planar flow. To these, we add 
the further assumptions: (v) the flux F is constant with 
z, but gravity g may depend arbitrarily on z. This can 
serve to model winds from thin, isothermal accretion 
disks, in which case g grows first linearly with z and at 
large distances drops off as Alternatively, with F 
and g being constants, one can model the launch region 
of a stellar wind, but a well-known degeneracy occurs 
here (Poe, Owocki, & Castor 1990). (vi) Zero sound 
speed is assumed, a — 0, and (vii) we fix a = 1/2. 
Note that the Sobolev line force is independent of vth 
and therefore of a. LDWs are hypersonic, and except 
near the photosphere, gas pressure plays no role. The 
1-D continuity and Euler equations are. 
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with constant Cq. We consider first stationary solu- 
tions, pv = const. A normalized quantity m = pv j p^Vc 
is introduced, where pjyz) and v^iz) are the density and 
velocity law of the critical wind, which is defined be- 
low. Besides m, a second, new hydrodynamic variable, 
w' — vv' , is defined, and the Euler equation becomes, 



+ g{z) ~ Cy/w'/m = 0. 



(6) 



The flux F was absorbed into the constant C. At each 
z, is a quadratic equation in \Aj7, with solutions, 

1 



w = 



C ± - 4gm 



(7) 



The velocity law v{z) is obtained from vu' by quadra- 
ture. Solutions with ' are termed shallow, those with 
'+' are termed steep. If m < 1 (see below), shallow 
and steep solutions exist from z = Otooo. Ifm>l, 
shallow and steep solutions become imaginary in a cer- 
tain z interval. In this region, the line force ^ C/-\/m 
cannot balance gravity g. 

Within the family of shallow wind solutions, the 
mass flux pv increases monotonically with terminal 
speed, while for steep solutions the trend is opposite. 
The largest mass flux which keeps the solution every- 
where real defines the critical wind, mc = 1. Setting 
the square root in eq. (^ to 0, this implies for C, 



(8) 



where gc = g{zc) means gravity at the critical point of 
the critical solution. How is Zc found? Differentiating 
the stationary Euler equation, E(^z,w'{z)'j = 0, with 
respect to z and using dE/dw'^ = at the critical point 
(crossing of solutions), one finds 
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dZr 



dg 
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(9) 



Hence, the critical point coincides with the gravity 
maximum. This is not an accident, but expresses that 
the critical point lies at the bottleneck of the flow, as 
for a Laval nozzle (Abbott 1980). If the flux F varies 
with z, the generalized area function depends also on 
F, and the critical point no longer coincides with the 
gravity maximum (see Feldmeier & Shlosman 1999 for 
examples). If 5 = const, the critical point degenerates, 
and every point in the flow becomes critical. In stellar 
wind calculations, the correct critical point location is 
found by including the finite cone correction factor for 
the stellar disk as an 'area' function (Pauldrach, Puis, 
& Kudritzki 1986; Friend & Abbott 1986). 
The wind solution becomes, 
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(10) 



At the critical point, shallow and steep solutions with 
m — 1 merge in such a way that the slope in passing 
from one to the other is continuous. Staying instead 
on either shallow or steep solutions introduces a dis- 
continuity in v"{zc). Discontinuities in derivatives of 
hydrodynamic variables, termed weak discontinuities, 
lie on flow characteristics (Courant & Hilbert 1968). 
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Characteristics are the space-time trajectories of wave 
phases. Indeed, we find below that the critical point is 
a barrier for Abbott waves. 

It is at this point that the question arises, which 
solution the wind adopts — a shallow, steep, or critical 
one. This issue will be resolved by discussing runaway 
of shallow solutions. We shall find that Abbott waves 
are the prime driver of this evolution. 

2.2. Green's function 

We derive the Green's function for Abbott waves 
in Sobolev approximation. The Green's function gives 
the response of a medium to a localized, delta function 
perturbation in space and time, and is complementary 
to the harmonic dispersion analysis of Abbott (1980) 
and Owocki & Rybicki (1984). Since localized pertur- 
bations contain many harmonics, a Green's function 
describes wave interference. This is clearly seen for 
surface water waves, whose Green's function is known 
from Fresnel diffraction in optics (Lamb 1932, p. 386). 
For simplicity, we consider only a single, optically thick 
line, with Sobolev force (a = 1), 
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dv 

dz 



(11) 



Density p was absorbed into the constant A. We 
assume WKB approximation to hold (slowly varying 
background flow) and consider velocity perturbations 
only. The characteristic analysis in the next section will 
show that the Abbott wave amplitude is v' / p, hence, 
Abbott waves are not annihilated by this restriction to 
velocity perturbations. The linearized Euler equation 
for small perturbations is 



d_ 

dt 



Sv{z,t) = 5gi{z,t) = A5v'{z,t). 



(12) 



The Green's function problem is posed by specifying as 
initial conditions, 

5v{z,Q) = 5{z- zq). (13) 
Multiplying eq. (|l2|) by e~''^^ and integrating over z. 



d_ 

dt 



Sv{k,t) = ikA6v{k,t), 



(14) 



where a bar indicates Fourier transforms, Sv = J Sv 
e~^^^-dz. The right hand side was obtained by inte- 
gration by parts, assuming (5w(— oo,t) — 6v{oo,t) ~ 0. 
This is shown a posteriori. The solution of ( p^ ) is 

^(fc,t) = 6e*'='^*, (15) 

with constant b. Fourier transforming (^3|) from z to k 
space. 



and 

S^{k,t) = e*'=(^*-^o). 
Fourier transforming back to z space. 



(17) 



1 r°° 

Sv{z, t) = — dk e^fc-e^'^f^*-"") = 5{z - z^ + At). 

Therefore, the initial delta function propagates without 
dispersion towards smaller z, at an Abbott speed —A. 
Furthermore, 5v = Q dX z = ±oo, as assumed. Since 
no wave dispersion occurs, the same Abbott speed A is 
also obtained by considering harmonic perturbations. 
Inserting 5v = Sv e*('^^~'^*) in (^2|) gives as phase and 
group speed. 



U! doj ^ 
k dk 



(19) 



The Green's function, G, is defined by {F an arbitrary 
function) , 



F{z,t)^ dz' G{z-~ z\t)F{z\0). (20) 

-/ — OO 

From (ll) and (|l|), 

G{z,t)^S{z + At), (21) 

a result first obtained by Owocki & Rybicki (1986). 

The present case of optically thick lines only corre- 
sponds to a = 1. An explicit expression for the Abbott 
speed is not relevant then: opposed to all cases a < 1, 
a — 1 poses no eigenvalue problem for to. We return 
therefore to a = 1/2. 

2.3. Abbott wave characteristics 

Besides a harmonic and Green's function analysis, a 
characteristic analysis can be given for Abbott waves. 
Especially, the latter is not restricted to linear waves. 
Inserting C from (||), the equations of motion (4, 5) 
become (dots indicate time derivatives). 



p + vp' + pv' ~ 




where we introduced the constant 



r = ^/gcPcVc 



(24) 



Sv(k,0) = e-''''" = b, 



(16) 



To bring these equations into characteristic form, we 
first write the continuity equation formally as K(p, v) = 
0, the Euler equation as E{p, v) = 0. For non-linear, 
first order systems of partial differential equations. 
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K = and E = 0, in two unknown variables p and 
V, the latter being functions of coordinates t and z, the 
characteristic directions or speeds, a, are determined 
by (Courant & Hilbert 1968; vol. II, chap. 5, §3) 



The Euler equation in characteristic form is therefore. 



^aKp + Kp 
-aEo + En' 



-aEy + Ey' 



= 0, 



(25) 



where Kpi 



dK/dp', Ey = dE/dv, etc. We use the 
symbol a, hitherto reserved for the sound speed, also 
for characteristic speeds. The meaning should be clear 
from the context. Inserting K and E in (|2E 



hence, 



-a + v 




P 



= 0, 



A 



/ pv' 



(26) 



(27) 



in the observers frame. The Abbott speed is again 
denoted A. In the comoving frame, a+ = and 
a_ = —V I \J pv' . The downstream (+), slow wave mode 
corresponds to sound waves. The upstream (— ), fast 
mode is of radiative origin. 

A simple, heuristic argument can be given for the 
occurence of Abbott waves. Consider a long scale per- 
turbation of a stationary velocity law. At the node 
where the velocity gradient gets steepened, the Sobolev 
line force increases. The gas is accelerated to larger 
speeds, hence the node shifts inwards. Similarly, the 
node where the velocity law becomes shallower shifts 
inwards. The node shift corresponds to phase propaga- 
tion of a harmonic wave. 

Next, we bring the equations of motion into char- 
acteristic form. To this end, the Euler equation is 
quasi-linearized by differentiating it with respect to z 
(Courant & Hilbert 1968), introducing a new, funda- 
mental variable f = v' , 

f + vf + f-^(f'-f^)+g'^0. 

Vfp\ Pj (28) 

Re-bracketing and multiplying with p. 



pf + pAf + pf + -^fp' + pg' = 0. 



(29) 



Using the continuity equation. 



= pf + pAf - fp-[v--j=],fp' + pg' 



^pf + pAf -fp- Afp' + pg' 
= p2 ^dt + Ad,)^+pg'. 



(30) 



{dt + Ad,)^ 



P 



1. 

P 



(31) 



with A from (P7|). We assume that WKB approxima- 
tion applies, i.e., that the temporal and spatial deriva- 
tives on the left hand side are individually much larger 
than the right hand side, hence the latter can be ne- 
glected. In a frame moving at speed —A, the func- 
tion v' / p is constant, and can be interpreted as a wave 
amplitude. Note that v' / p is inversely proportional to 
the Sobolev line optical depth, indicating that Abbott 
waves are indeed a radiative mode. 

Introducing / in the continuity equation puts it into 
characteristic form. 



{dt + vd^)p= - fp. 



(32) 



Here, fp is an inhomogcneous term. WKB approxi- 
mation cannot be assumed here, since / may vary on 
short scales. The wave amplitude, p, is no longer con- 
stant along V characteristics, but changes according to 
this (ordinary) differential equation. Since gas pres- 
sure p scales with density, this equation shows that the 
outward mode corresponds to sound. 

For stationary winds, the Abbott speed in the ob- 
servers frame becomes {a — 1/2), 



A^vil- 



mw 



= u 1 - 



1 



^mw'/w'J ' (33) 



since w'^ — gc from (|, |). With A^ = 0, the criti- 
cal point is a stagnation point for Abbott waves. For 
shallow winds, m < 1 and w' /w'^ < 1, hence ^ < 
and Abbott waves propagate upstream from any z to 
the photosphere located at z = 0. Shallow solutions 
are, therefore, sub-abbottic, and are the analog to so- 
lar wind breezes. For steep solutions, y^mw' /w'^ > 1 
from (p^. Hence A > 0, and steep solutions are super- 
abbottic. Once the wind has adopted a steep solu- 
tion, the flow can no longer communicate with the wind 
base: steep solutions cannot evolve by means of Abbott 
waves. 

2.4. Negative velocity gradients 

A surprising result occurs when we allow for negative 
velocity gradients, v' < 0, somewhere in the wind. This 
corresponds to flow deceleration, not necessarily to ac- 
cretion instead of wind. The Sobolev force is blind to 
the sign of the velocity gradient. All that is important 
is the presence of a velocity gradient, to Doppler-shift 
ions out of the absorption shadow of intervening ions. 
Hence, a natural generalization of the Sobolev line force 
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is, 

5i = 2rW^. (34) 

V P 

This liolds for a purely local force. However, if u' < 0, 
the velocity law is non-monotonic and multiple reso- 
nance locations occur. Radiative transfer is no longer 
local, as photons are absorbed and scattered at dif- 
ferent locations. The incident flux is no longer de- 
termined by the photospheric flux F alone, but for- 
ward and backward scattering has to be accounted for. 
The constant C becomes frequency and velocity depen- 
dent. Rybicki & Hummer (1978) introduced a general- 
ized Sobolev method for non-monotonic velocity laws, 
where the radiation field is found by iteration. This in- 
troduces interesting, non-local effects into Abbott wave 
propagation (action at a distance). We postpone such 
an analysis to a future paper, and proceed here in a 
simpler fashion. Together with gi ^ the oppo- 

site extreme gi ^ ■\/max(u', 0) is treated. In the latter 
force, all radiation is assumed to be absorbed at the 
first resonance location, where necessarily v' > 0. The 
line force according to Rybicki & Hummer (1978) lies 
in between these two extremes. 

Repeating the above steps for these generalized line 
forces, the Euler equation maintains its characteristic 
form, 

(dt + Ad,)- = -^, 
P P 

with Abbott speed A = v i — for v' > 0, and 

V pv' 

iv+ / , if.gi-/H 

[v if 51 - v/max(w',0) (35) 

for v' < 0. Therefore, if the velocity gradient is nega- 
tive, Abbott waves propagate downstream, with a posi- 
tive (or zero) comoving frame velocity along all solution 
types, whether shallow, steep, or critical. This is pe- 
culiar, since Abbott waves appeared so far as upstream 
mode. (Note that, for gi ^ max(w', 0), the line force 
drops out of the Euler equation if v' < 0. Both wave 
modes become ordinary sound then.) We conclude that 
regions with v' < cannot communicate with the wind 
base. 

3. Abbott wave runaway 
3.1. Method 

In the remainders, we study wave propagation in 
LDWs numerically, using a standard time-explicit, Eu- 
lerian grid code (van Leer advection on staggered grids) . 



Non-refiecting Riemann boundary conditions for Ab- 
bott waves are used (Feldmeier et al. 2001). As inner 
boundary, z = 0.1 is chosen to avoid negative speeds 
when Abbott waves leave the mesh at the wind base: 
numerical artefacts may result when the v character- 
istic changes its direction. For gravity, we assume 
g — z/{l + z'^), with a maximum at z = 1 . Since g and z 
are normalized, so are v and t. Steep solutions are of no 
further interest here, since they are super-abbottic and, 
therefore, numerically stable. Furthermore, since they 
start supersonically at the wind base, they are unphys- 
ical. We are left with shallow winds, which can evolve 
towards the critical solution by means of Abbott waves. 
Shallow solutions are numerically unstable if pure out- 
flow boundary conditions are used. The mechanism of 
this runaway is not easy to analyse, because of numeri- 
cal complications in the vicinity of the boundary, where 
the nature of the difference scheme changes. 

To clearly separate effects of boundary conditions 
from wave dynamics, we introduce controlled, explicit 
perturbations in the middle of the calculational do- 
main. 

3.2. Mechanism of the runaway 

Figure |] demonstrates the result derived above, that 
positive and negative velocity slopes propagate in op- 
posite directions. The initial conditions are a shallow 
velocity law which is perturbed by a triangular wave 
train. The subsequent evolution of this wave train fol- 
lows from the kinematics of the velocity slopes. In strict 
mathematical terms, the plateaus which form in Fig. |l| 
correspond to centered rarefaction waves. We postpone 
such an analysis to a forthcoming paper. The essential 
result from the flgure is that the wind speed in the saw- 
tooth evolves asymmetrically, towards larger values. 

The whole sawtooth pattern moves upstream as an 
Abbott wave. This is of a prime importance for our 
understanding of the observed runaway. 

Namely, if a perturbation is fed into the wind con- 
tinuously over time, the whole inner wind is eventually 
lifted towards larger speeds and mass loss rates. The 
same is true for the outer wind, first directly by the run- 
away, and second as a consequence of the accelerated 
inner gas propagating outwards. 

We consider a coherent, sinusoidal velocity pertur- 
bation of period P and maximum amplitude S, which 
is fed into the flow at a flxed location z. The funda- 
mental hydrodynamic variables used in the code are p 
and pv. After each time step St, perturbations 

pv pv + pSv, p p{l — Sv/\A\), 

(36) 
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Fig. 1. — Evolution of a triangular wave train in a shal- 
low wind (the latter shown as dotted line). For clarity, 
the velocity law is plotted with a negative, vertical off- 
set which increases with time. The dashed lines show 
the shallow solution with m = 1 and the critical solu- 
tion. 




Fig. 2. — Left panel: stable Abbott wave propagation 
along a shallow velocity law. A sinusoidal perturbation 
with amplitude S = 0.04 and period P = 1 is applied 
at Zo = 2. Right panel: Abbott wave runaway if the 
amplitude is doubled to S" = 0.08. 

with 

, 2ttS 2Trt , , 

ov — St cos (37) 

are applied to pv and p on a single mesh point. The 
density fluctuations follow from the continuity equa- 
tion, 5p/ p w — (5u/|A|. For linear waves, the observers 
frame Abbott speed is A « —1.05 at = 2 for m = 0.8. 

For sufficiently small amplitudes S, v' remains pos- 
itive, and Abbott waves propagate in a stable fashion 
towards the wind base. This is shown in the left panel 
of Figure ||, for P = 1 and S = 0.04. Doubling the per- 



turbation amplitude to S' = 0.08 implies wind runaway 
towards the critical solution, as is shown in the right 
panel of Figure |^. The wind converges to m = 1 ev- 
erywhere (not shown). Instead of adopting the critical, 
accelerating branch, the velocity law jumps at z > 
to the decelerating branch. This is even true for the 
converged, stationary solution as t 00. (The veloc- 
ity slope is so mildly negative for t —^ 00 that the wind 
speed is almost constant above — 2.) We add some 
further remarks on this issue below. 

The runaway results from the occurence of nega- 
tive velocity gradients. During excitation phases where 
w' < 0, the resulting line force perturbations are not 
sufficiently negative to compensate for positive line 
force perturbations during phases where v' > 0. Net 
acceleration of the wind results over a full excitation 
cycle is the consequence. 
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Fig. 3. — The runaway time series from Fig. y, shown 
with a constant vertical displacement between snap- 
shots. Dashed lines show phases where negative ve- 
locity perturbations are applied at Zo, leading to self- 
annihilating slopes. 

Figure ^ shows again the runaway time series of 
Fig. |[ with subsequent snapshots displaced vertically 
for clarity. During the negative perturbation half-cycle, 
— (5w(zo), negative velocity slopes propagate outwards 
from below Zg) and merge with positive slopes propa- 
gating inwards from above Zq. The merging slopes mu- 
tually annihilate. After a time P/2, the velocity law is 
left largely unaltered in presence of the perturbation. 
This causes the dense spacing of curves in Figure ^, 
especially at z > Ze, once every perturbation cycle. 
On this rather flat velocity law, a positive perturba- 



7 



tion +Sv is added during the next half-cycle. Here, 
inner, positive slopes propagate inwards and separate 
from outer, negative slopes which propagate outwards. 
Obviously, runaway is caused by positive perturbations, 
but its deeper origin is that negative perturbations are 
self-annihilating, and cannot balance positive ones. 

The runaway terminates when comes to lie on 
the critical solution, and Abbott waves can no longer 
propagate inwards. 

The velocity gradient at z > Zo is then still nega- 
tive, and in most of our simulations remains negative 
at all times: at Zo, the wind jumps to the decelerating 
branch with m — 1, which causes a kink in v{z). Run- 
away perturbations above Zc would be required to es- 
tablish a critical, accelerating velocity law in the outer 
wind. For certain combinations of model parameters, 
we find instead a critical solution with v' > over the 
whole mesh. The near-plateau above z^ evolves then to- 
wards larger speeds, and the velocity kink propagates 
outwards, eventually leaving the mesh. We believe that 
this is an artefact caused by boundary-induced, numer- 
ical runaway (Feldmeier & Shosman 2001). The latter 
can even occur when non-reflecting Abbott boundary 
conditions are used, via non-WKB (standing?) waves; 
our Riemann boundary conditions were formulated to 
annihilate WKB waves only. We prefer the situation 
w' < for z > Zo over any numerical runaway which 
would assist the present, physical runaway in reaching 
the critical solution. A future analysis of outer bound- 
ary conditions has to clarify this issue of the outer wind 
velocity law. 

We can still test the stability of the full, critical so- 
lution, as is done in Figure ^. Above the critical point, 
perturbation phases with v' > Q and v' < combine 
here, and propagate outwards as a smooth, marginally 
stable Abbott wave. 

3.3. Stationary overloaded winds 

So far, we assumed that wind perturbations are lo- 
cated above the critical point. We consider now the 
opposite case, Zc < z^. 

Figure ^ shows the wind velocity law resulting from 
a perturbation location Zo = 0.8 on the critical solution. 
A period P = 0.3 was chosen. A stationary wind with 
a broad deceleration region, v' < 0, develops. Abbott 
waves propagate outwards from Zo through the deceler- 
ating wind. Since the CAK critical point, which is the 
bottleneck of the flow, lies in the deceleration regime, 
the present solution should be overloaded. Indeed, the 
mass loss rate is found to be m = 1.05. 

This result is readily understood. The flow at Zc 
is still sub-abbottic, because Zc < z^. Runaway to 
larger v occurs until inward propagation from Ze be- 
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Fig. 4. — The critical, CAK solution is stable with 
respect to Abbott waves excited above the critical 
point. Perturbation phases with positive and negative 
v' combine to an outward propagating, harmonic Ab- 
bott wave. 

comes impossible. The velocity at z < Ze is then ev- 
erywhere larger than the critical speed. The runaway 
stops when the inward Abbott speed becomes zero in 
the observers frame. We are interested in a stationary 
solution, hence, A = in (p^), implying mw' — g^- The 
square root in ( |lo|) has to vanish then, hence m — gd g. 

So far, we restricted ourselves Xo g ~ and m — \. 
But critical points can also occur along overloaded so- 
lutions. In this case, m = Qc/ g{zc) > 1, and the square 
root in ( ]l0| ) becomes imaginary at z > Zq. Remem- 
ber that gc refers to gravity at the critical point of the 
critical solution, which is the maximum of g. Hence, 
g{zc) < gc for any other possible critical point. The 
perturbation site Zc stops communicating with the wind 
base only when it becomes a critical point itself. At Zo, 
the overloaded wind jumps to the decelerating branch, 
resulting in a kink in the velocity law. This kink has 
all the attributes of a critical point. Especially, as an 
Abbott wave barrier, it shuts off communication with 
the wind base. 

It is easily seen that a critical point at z < Ze would 
imply that Zc is already super-abbottic, which cannot 
happen during runaway. Since g{z) is a monotonically 
growing function with a maximum gc, the perturba- 
tion site Zc < Zc is the first one that becomes critical 
to Abbott waves on the evolving overloaded solution. 
We emphasize that the amount of overloading is solely 
determined by Zg. 
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1999). Another jump to the decelerating branch is un- 
avoidable. Since the local wind speed is much larger 
than the local escape speed, the velocity law stays es- 
sentially flat above this kink. 



3.4. Time-dependent overloaded winds 

Already a minor increase in the mass loss rate causes 
a broad deceleration region. This is a consequence 
of g{z) having a broad maximum at Zc- Integrating 
w'{z) = vv' numerically, we find that v{l/zc) = for 
Zc = 0.660. If Zc is still smaller, the gas starts to fall 
back towards the photosphere before it reaches l/ze, 
and collides with upwards streaming gas. A station- 
ary solution is no longer possible. Figure ^ shows, 
for Ze = 0.5, that an outward propagating sequence 
of shocks and dense shells forms. The shock spacing is 
determined by the perturbation period. 



Fig. 5. — Velocity law for a wind with harmonic per- 
turbation at Ze = 0.8, below the critical point Zc = 1. 
The wind converges to a stationary, overloaded solu- 
tion. The dotted line shows the initial conditions, the 
critical CAK solution. 



At some height above the gravity maximum, the 
decelerating solution jumps back to the accelerating 
branch, giving a second kink: the wind has overcome 
the bottleneck, and starts to accelerate again. In the 
present wind model, the (imaginary, overloaded) solu- 
tion becomes real again at z = 1/zo (zo being the loca- 
tion where the solution becomes first imaginary). For 
Ze = 0.8, the overloading is m = 1.025 and the wind 
can start accelerating again at z = 1.25. In the simu- 
lation however, m = 1.05, and the second kink occurs 
at z = 1.6. The discrepance in m can be attributed 
to mesh resolution, which blurrs out z^. This leaves 
still unanswered why the wind 'waits' so long before it 
starts accelerating again. 

Furthermore, starting with a shallow wind as initial 
conditions instead of a critical wind, the simulation 
converges again to an overloaded solution. However, 
the jump from the decelerating back to the accelerat- 
ing branch does not terminate on the steep, but already 
on the shallow overloaded solution with m = 1.05. This 
is a numerical artefact caused by non-reflecting, outer 
boundary conditions, which try to maintain shallow so- 
lutions. As above, we prefer this situation over any 
boundary-induced, numerical runaway. 

A third kink may occur. For realistic radiative fluxes 
above accretion disks, we find that at still larger dis- 
tances from the disk, the line force drops below absolute 
gravity for all solution types (Feldmeier & Shlosman 
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Fig. 6. — Velocity v and mass flux to of a wind per- 
turbed at Zo = 0.5. During each perturbation cycle, gas 
falls back towards the photosphere and collides with 
upwards streaming gas. Shocks form, and propagate 
outwards. 

As a technical comment, we add that the latter sim- 
ulation tends to develop extremely strong rarefactions 
and, correspondingly, very large velocities. The latter 
cause the Courant time step to approach zero. This is a 
well-known artefact of the power law form of the CAK 
line force, gi (x (w'/p)", implying ^ oo as p ^ 0. 
The true line force reaches instead a finite maximum 
in a flow which is optically thin even in very strong 
lines (Abbott 1982). This is achieved by truncating the 
line distribution function (Owocki, Castor, & Rybicki 
1988). Using a simpler approach which suffices for the 
present pvirposes, we simply truncate m at large values. 
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Even for — 0.660 implying negative speeds, over- 
loading is small, only 9% above the critical CAK value 
(m = 1.088). This shows that the CAK mass loss rate is 
a significant upper limit for mass loss rates in LDWs, 
with or without mass overloading. If overloading oc- 
curs, it should be detectable only via broad, decelerat- 
ing flow regions. 

It is important to distinguish the present runaway 
from the well-known line-driven instability (Lucy & 
Solomon 1970; Owocki & Rybicki 1984; Lucy 1984), 
which leads to strong shocks and dense, narrow shells 
in the wind. For perturbations longer than the Sobolev 
scale, the instability can be understood as a linear pro- 
cess in second-order Sobolev approximation, including 
curvature terms v" . By contrast, runaway occurs al- 
ready in first-order Sobolev approximation, yet, finite 
perturbations are required to achieve v' < 0. 

3.5. Triggering perturbations 

Amplitudes. Still, even perturbations with small am- 
plitudes can lead to runaway, if their wavelength is 
sufficiently short to cause negative v'. One expects, 
however, that dissipative processes like heat conduction 
prevent runaway for short wavelength perturbations. 

Coherence. In a hydrodynamic instability, a feed- 
back cycle amplifies small initial perturbations. Even 
perturbations of short duration trigger unstable growth. 
By contrast, the present wind runaway requires a con- 
tinuously maintained perturbation seed. No feedback 
occurs. The runaway is a pure wave feature caused by 
the peculiar w'-asymmetry of the line force. 

Perturbations of short duration lead to a localized 
runaway. The wind is shifted towards larger v and m 
within a small z interval. If the perturbation ceases, 
so does the runaway. The region of increased v and 
TO propagates upstream as an Abbott wave, and even- 
tually lifts the wind base to a stable shallow solution 
with slightly higher m. Successive Abbott waves lift 
the wind until the critical or an overloaded solution is 
reached. 

Surrounding medium. One expects that mismatches 
at the outer boundary, where the wind propagates 
into a medium of given properties, create perturbations 
which can trigger runaway. This is supported by the 
fact that at large z the velocity law is almost flat, and 
even small-amplitude perturbations can cause v' < 0. 
The Abbott speed is large at large z, since A ~ 1/Vi7, 
hence perturbations propagate inwards quickly. We 
performed tests with the outer boundary at z = 40. 
Inherent interpolation errors are sufficient then to set 
off boundary runaway. 

Other sources. Inner wind perturbations could oc- 
cur due to prevalent shocks from the line driven in- 



stability (Owocki et al. 1988). The critical point lies 
usually close to the sonic point (Owocki & Puis 1999; 
Feldmeier & Shlosman 1999), and it is therefore not 
easy to contemplate strong shocks below the critical 
point. This means that perturbations at large and not 
at small z dominate the runaway, and drive the wind 
to a critical solution. Hence, overloaded winds formed 
by internally-generated perturbations should be rare. 

3.6. Critical points and mass loss rate 

We discuss here some general issues related to criti- 
cal points and mass loss rates in different types of stellar 
winds. 

Holzer's wind laws. The idea that upstream, in- 
wards propagating waves adapt the wind base to outer 
(boundary) conditions is fundamental to solar wind 
theory. From this grew the recognition that an out- 
ward force which is applied above the sonic (critical) 
point does not affect the mass loss rate, but only ac- 
celerates the fiow; whereas a force applied below the 
sonic point increases the mass loss rate, but has a van- 
ishing effect on terminal wind speeds (Leer & Holzer 
1980; Holzer 1987). These 'wind laws' have proven to 
be of interest for cases far beyond the coronal winds for 
which they were first applied. We refer to Lamers & 
Cassinelli (1999) for a detailed discussion. 

Most relevant to us is the case of a force applied 
in the subcritical wind region, causing enhanced mass 
loss. For coronal winds, the subsonic region has es- 
sentially a barometric density stratification. Any extra 
force assisting the pressure gradient helps to establish 
a larger scale height. From the continuity equation, a 
shallower velocity gradient results, though the termi- 
nal speed is hardly affected. The corresponding situ- 
ation, analyzed by us for LDWs, is even simpler. No 
outward force occurs besides line driving. Our choice 
of zero sound speed emphasizes that the barometric 
stratification plays no role for mass loss runaway or for 
establishing an overloaded solution. Overloaded LDW 
solutions have a steeper velocity gradient than the crit- 
ical solution, as can be seen from eq. (|7|). Physically, a 
steeper velocity gradient is required to create a gener- 
alized critical point below the CAK critical point. In 
LDWs, a critical acceleration w' = vv' , not a critical 
speed V is adopted at the critical point, and prevents 
further Abbott wave runaway. 

Physical relevance of critical points. The physical 
relevance of the CAK critical point was questioned by 
Lucy (1975, 1998), arguing that it may be an arti- 
fact of Sobolev approximation. Holzer (1987, p. 296) 
doubts that critical points, i.e., singularities of the flow 
equations, generally coincide with those points beyond 
which "information relevant to the acceleration of the 
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wind" can no longer be transported upstream. We have 
tried in the present paper to re-estabhsh a more tradi- 
tional viewpoint (Courant & Friedrichs 1948; Courant 
& Hilbert 1968), namely that a critical point as a math- 
ematical singularity leaves certain derivatives of flow 
quantities undetermined. This is only possible along 
characteristics, which are the space-time trajectories of 
the specific waves in the problem. This establishes the 
transonic (or trans-abbottic, etc) nature of the critical 
point. Similarly, one can argue: the critical point is 
a singularity because different solution branches cross 
there. Hence, by passing continuously from one branch 
to another (actually: by staying on a certain branch), 
higher-order or weak discontinuities appear in the flow 
properties; for LDWs, in v" . Again, weak discontinu- 
ities propagate at characteristic wave speed, making 
the critical point transonic. 

Mass loss rate as an eigenvalue. We remind the 
reader of a deep difference between coronal and line- 
driven winds. For coronal winds, the mass loss rate is a 
free parameter within wide margins, and is determined 
by the base density. For LDWs, on the other hand, the 
mass loss rate is a unique, discrete eigenvalue (CAK). 
This is a consequence of the line acceleration depending 
on p, see Lamers & Cassinelli (1999). 

Abbott speed and speed of light. In the literature, 
the propagation speed of radiative waves is occasion- 
ally identified with the speed of light. This is not true 
in general: the Abbott speed and, in magnetized fiows, 
the Alfven speed are smaller than the speed of light. 
To make the point totally obvious, note that the diffu- 
sion speed of photons through stellar interiors is much 
smaller than the speed of light. The basic cause is here 
the same as for Abbott waves: optical depths larger 
than unity. 

3.7. Deep-seated X-rays, infall, and mass-over- 
loading 

Besides adapting the wind base to outer flow condi- 
tions, wc find that Abbott waves can even lead to mass 
infiow. This seems like a novel feature of these waves 
and it leads to interesting observational consequences, 
like explaining the formation of hard X-rays anoma- 
lously close to the surfaces of hot stars. Formation 
depths of r « 1.1 i?* are deduced from X-ray emission 
lines observed with the Chandra satellite (for C Ori: 
Waldron & Cassinelli 2000; for ( Pup: Cassinelli et 
al. 2001). The favored model for X-ray emission from 
hot stars is via strong shocks (Lucy 1982) from the de- 
shadowing instability (Lucy & Solomon 1970; Owocki 
& Rybicki 1984). The shocks become strongest in col- 
lisions of fast clouds with dense gas shells (Feldmeier, 
Puis, & Pauldrach 1997). From theoretical arguments 



(line drag effect: Lucy 1984) and numerical simula- 
tions it appears that this shock scenario cannot explain 
X-rays originating from the small heights mentioned 
above. 

Howk et al. (2000) suggest gas infall as possible ori- 
gin of X-rays from near the photosphere. They consider 
a ballistic model of stalled wind clouds falling back to- 
wards the star, to explain the hard X-rays observed 
from r Sco (Cassinelli et al. 1994). At any instant, 
~ 10% of the mass flux from the star resides in nu- 
merous, > 1000, clouds which barely reach distances 
of r = 2i?* before falling back inwards. The overall 
resemblance of this cloud model with time-dependent 
wind overloading as shown in Fig. ^ is striking. In the 
latter case, the overloading is sufficiently mild that up- 
streaming gas can push downfalling clumps outwards 
(cf. the drag force of ambient wind gas in Howk et 
al.). Stronger overloading is achieved by shifting the 
Abbott wave source further inwards, but cannot be ad- 
dressed using our simple numerical approach. Clearly, 
two-dimensional simulations are required to model true 
infall. A systematic study relating the Howk et al. ap- 
proach with ours is in preparation. 

Note that the critical point in stellar wind models 
accounting for the finite cone effect lies at r < 1.1 
(Pauldrach et al. 1986). An overloaded wind starts 
to decelerate below the critical point, and may al- 
ready reach negative or infall speeds at similarly small 
heights. This shows the relevance of overloading in 
understanding the Chandra observations mentioned 
above. The topic of stalling and backfalling gas has 
aquired some general attention in recent papers on 
LDWs, and is discussed in different contexts by Friend 
& Abbott (1986), Poe et al. (1990), Koninx (1992), 
Proga, Stone, & Drew (1998), and Porter & Skouza 
(1999). 

We finally mention another idea related to over- 
loaded flows, which has played some role in sharpen- 
ing our understanding of coronal winds. Cannon & 
Thomas (1977) and Thomas (1982) challenged Parker's 
(1958, 1960) theory of the solar wind. They postu- 
lated an outwards directed, sub-photospheric mass fiux. 
In analogy with a Laval nozzle (e.g.. Chapman 2000, 
Chapter 7) into which a too strong mass or energy in- 
flux is fed (Cannon & Thomas 1977 use the terminology 
'imperfect nozzle' instead of 'overloading'), the solar 
outflow should choke, creating shocks below the nozzle 
throat, or the critical point of the smooth (perfect noz- 
zle) flow. The shocks are responsible for heating the 
chromosphere and corona, making the latter a conse- 
quence, not the origin of outflow from Sun. This model 
was ruled out by Parker (1981) and Wolfson & Holzer 
(1982). 
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Still, it leads over to another interesting difference 
between coronal and Laval nozzle flows on the one side 
and LDWs on the other: the former two do not al- 
low for continuous, overloaded solutions; instead, over- 
loaded flows have shocks in the vr plane where the crit- 
ical point resides. By contrast, overloaded LDWs show 
jumps in the critical vv'-r plane, i.e., jumps in the wind 
acceleration which correspond to kinks in a continuous 
velocity law. 

4. Summary 

We have studied the stability of shallow and steep so- 
lutions for accretion disk and stellar line driven winds. 
This was done by introducing flow perturbations into 
the wind at a fixed height. For sufficiently large per- 
turbation amplitudes, negative velocity gradients occur 
in the wind and cause runaway towards the critical so- 
lution. The origin of this runaway is an asymmetry in 
the line force: negative velocity gradients cause a force 
decrease which cannot balance the force increase dur- 
ing phases where the velocity law gets steepened. Net 
acceleration results over a full perturbation cycle. 

A new type of waves, termed Abbott waves, are ex- 
ited by the perturbations. They provide a communi- 
cation channel between different parts of the wind and 
define an additional critical point in the flow, down- 
stream from the sonic point. Shallow solutions are 
subcritical everywhere, steep solutions are supercriti- 
cal. Along shallow solutions, Abbott waves propagate 
upstream towards the photosphere at a high speed as 
a radiative mode, and creep downstream at a fraction 
of the sonic speed as an acoustic mode. The critical 
solution is the one that switches continuously from the 
shallow to the steep branch at the critical point. 

Inward propagating Abbott waves turn the local, 
asymmetric response to velocity perturbations into a 
global runaway towards the critical solution. The con- 
verged, steady wind solution depends on the location of 
seed perturbations. Three spatial domains can be dis- 
tinguished, (i) If perturbations are prevalent in outer 
wind regions, above the critical point, the wind set- 
tles on the critical solution, (ii) If perturbations occur 
just below the critical point, runaway proceeds to a 
stationary, so-called overloaded solution. Such a flow 
is characterized by a broad deceleration domain. The 
line force cannot balance gravity in a vicinity of the 
critical point (the bottleneck of the flow) because of a 
supercritical mass loss rate, (iii) Perturbations close to 
the photosphere result in the wind being decelerated 
to negative speeds. Gas falls back towards the pho- 
tosphere, and collides with upwards streaming gas. A 
time-dependent, overloaded wind results with a train 
of shocks and shells propagating outwards. 



The runaway nKxiianisni discussed here depends 
solely on the asymmetry of the line force with respect 
to velocity perturbations which cause local flow decel- 
eration, dv/dr < 0. It should, therefore, be a robust 
feature and not depend on Sobolev approximation. 
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